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We investigate numerically the magnetisation dynamics of an array of nano-disks interacting through the 
magneto-dipolar coupling. In the presence of a temperature gradient, the chain reaches a non-equilibrium steady 
state where energy and magnetisation currents propagate. This effect can be described as the flow of energy 
and particle currents in an off-equilibrium discrete nonlinear Schrodinger (DNLS) equation. This model makes 
transparent the transport properties of the system and allows for a precise definition of temperature and chemical 
potential for a precessing spin. The present study proposes a novel setup for the spin-Seebeck effect, and shows 
that its qualitative features can be captured by a general oscillator-chain model 
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I. INTRODUCTION 

Off-equilibrium dynamics of classical and quantum many- 
particle systems is is a wide topic, ranging from the theoret¬ 
ical foundations to the development of innovative ideas for 
nanoscale thermal management. The nanoscale control of 
heat flows offer promising opportunities for novel energy har¬ 
vesting devices, and constitutes a fertile terrain to with possi¬ 
ble future applications to nanotechnologies and effective en¬ 
ergetic resources. 

In this context, simple models of classical nonlinear oscil¬ 
lators have been investigated to gain a deeper understanding 
of heat transfer processes far from thermal equilibrium [1- 
3]. In the case of systems admitting two (or more) conserved 
quantities, coupled transport (for instance of energy and mass) 
is an issue of basic relevance in connection with thermoelec¬ 
tric phenomena [4] whereby temperature gradients can be em¬ 
ployed to generate electric currents. 

A related topic concerns coupled transport in magnetic sys¬ 
tems. The recent discovery of the spin-Seebeck effect has 
opened the new field of spin-caloritronics [5-7], where the 
transport properties of thermally driven spin systems is in fo¬ 
cus. The propagation of spin wave (SW) current driven by 
a thermal gradient has been extensively studied in systems 
of exchange coupled spins, using the well established micro- 
magnetic formalism [8-11]. In the context of statistical me¬ 
chanics, although several studies of heat transport in classical 
Heisenberg spin chains exists [12, 13], this topic has not been 
treated in detail up to now. 

Here, we investigate transport in a novel setup, which con¬ 
sists of an array of magnetic nano-disks coupled through the 
magneto-dipolar interaction, and interacting with different ex¬ 
ternal reservoirs. As a result, energy and magnetisation cur¬ 
rents are carried by damped dipolar spin waves. 

To obtain a better insight, we follow a simplified physical 
picture: one can think of a spin system as a chain of nonlin¬ 
ear oscillators. Intuitively, a thermal bath acts as a random 


force on some part of the chain, whose effect is the propaga¬ 
tion of oscillations (spin waves) in the system, in the form of 
energy and magnetisation currents. More specifically, we will 
model the system as an open discrete nonlinear Schrodinger 
(DNLS) system that steadily exchanges energy with external 
reservoirs, a setup that has been considered only recently [14- 
lb]. In particular, we will employ a Langevin thermostatting 
scheme whereby complex white noises and dissipative cou¬ 
plings are added to the equation [17, 18]. Altogether, this 
guarantees that the chain reaches thermal equilibrium when 
put in contact with thermostats at the same temperature. As 
we will show, this approach has several advantages in terms 
of simplicity of description and it makes transparent the trans¬ 
port/thermodynamical properties of the micromagnetic sys¬ 
tem. 

This paper is organised as follows. In Section II we review 
the equation of motion of the magnetisation for a macrospin 
interacting with a reservoir. In Section III we outline the 
derivation of the simplified oscillator model and describe its 
relations with the micromagnetic equations. In Section IV we 
study the steady-state properties of a spin-chain made of ten 
disks coupled via dipolar interaction through micromagnetic 
simulations. Using the formalism of the DNLS equation, the 
currents and the effective spin temperature of the system are 
calculated. The results are compared with the direct simula¬ 
tion of the oscillator model in Section V. Finally, the conclu¬ 
sions of this work are contained in the last Section. 


II. PHYSICAL SYSTEM 

The system studied here, shown in Fig. 1 consists of an array 
of N = 10 identical nano-disks made of Permalloy (Py), cou¬ 
pled through the dipolar interaction. The first and last disks 
are coupled to Langevin thermal baths with temperatures re¬ 
spectively T+. Thermal fluctuations excite the SW modes 
of the system. In the presence of a temperature difference 
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Figure 1: Chain of disks coupled via dipolar interaction. Each disk 
behaves as a macrospin with precession frequency co n . The first and 
last element of each chain are coupled to thermal reservoirs with tem¬ 
peratures T +, that control the propagation of energy and magnetisa¬ 
tion currents. 


AT = T+ -T-, the disk chain reaches a non-equilibrium sta¬ 
tionary state where two coupled currents, of energy and mag¬ 
netisation, flow from the hot ( T +) to the cold (T_) reservoirs. 

The dynamics is investigated by means of micromagnetic 
simulations. The cornerstone of micromagnetism is the 
Landau-Lifschitz-Gilbert (LLG) equation of motion for a fer- 
romagnet [19-21] 


M = yM x H eff + —M x M, (1) 

M s 

which describes the precession of the local magnetisation vec¬ 
tor M(r,t ) around the effective field H c $. The first term of 
Eq.(l), proportional to the gyromagnetic ratio y, accounts for 
the precession. The second term describes energy dissipation 
at a rate proportional to the dimensionless Gilbert damping 
parameter a. In the absence of an external driving field, such 
as thermal fluctuations, spin transfer torque or rf fields, the 
magnetisation eventually aligns with iT e ff [21, 22]. The satu¬ 
ration magnetisation M s is the norm of the magnetisation vec¬ 
tor, which depends on the material and the geometry of the 
sample. The thermodynamical properties of the system are 
described by the Gibbs free energy 

T =ju 0 J' [AaM(t\ t) - M(r, t ) • H ext 

- l -M(r,t)-H &p ]d?r, (2) 

p o being the vacuum magnetic permeability. In the present 
case, Eq.(2) contains contributions respectively from ex¬ 
change energy, Zeeman interaction and dipolar interactions. 
The effective field in Eq.(l) is given by the functional deriva¬ 
tive = -(1 /yUo)#F/ 5M of the Gibbs free energy with 
respect to the magnetisation. It is the sum of the following 
three terms: 



Figure 2: Macrospin M precessing around the effective field H eff , 
described in terms of the stereographic projection if/. The phase 0 
describes the dynamics in the x-y plane, while the SW power p = \if /\ 2 
is related to the polar angle 6 (see text). 


along z. The exchange field is the short range interaction re¬ 
sponsible for the coherent precession of the magnetisation in¬ 
side each disk, while the applied field defines the precession 
axis. The dipolar stray field Eq.(3c) contains contributions 
from the volume charges pM = V • M(r', t) and the surface 
charges cr M = M(r' ,t)-n, where fi denotes the normal to the 
surface of the sample at point r'. The dipolar field acts as a 
demagnetising field in each disk and adds a coupling between 
the disks. It is responsible for the nonlinearity of the LLG 
equation, Eq. (1). 

Thermal fluctuations are introduced by adding to the effec¬ 
tive field the stochastic term 


H lh (r,t)= VDT(ri x ,7] y ,7] z ), (4) 

where rjj(r , t), j = (x, y, z), is a Gaussian random process with 
zero average and correlation (rij(r, t)Tjy(r\f^ = 8jyd(r - 
r')S(t - t'). [23, 24]. The term T = T(r) denotes the lo¬ 
cal temperature of the underlying phonon bath, while the cou¬ 
pling strength with the bath is 


D = 


2 ak B 

7P()VmM s 


( 5 ) 


Here k B is the Boltzmann constant and V M is the volume as¬ 
sociated to the magnetic moment M(r). 

In our micromagnetics simulations, where the sample is 
represented by a finite element tetrahedral mesh, the LLG 
equation, Eq. (1), is solved numerically at each mesh node. 
The coordinate r is discretised and corresponds to the posi¬ 
tions of the nodes, while V M corresponds the volume of each 
mesh elements. 


H exc =-AVM(r,t), (3a) 

H ex t =Hz, (3b) 


-f^dip — 


J_ r Pu(r-r') d3/i 1 r <r M (r-r') dV 
4 7r Jv \r - r'\ 3 4 7r J s |r-r'| 3 

(3c) 


III. COUPLED OSCILLATOR MODEL 

The dynamics of the chain is conveniently described by the 
volume-averaged magnetisation inside the ftth disk: 


Eqs. (3a) and (3b) are respectively the exchange field with 
exchange constant A and the external field with intensity H 


M n (t) = 



M(r n ,t)d 3 r n . 


( 6 ) 
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Due to the uniform precession of the magnetisation in each 
disk, the system can be modeled as an assembly of coupled 
macrospins M n (t). With this approximation, Eq.(l) can be 
written in the form of an equation of motion for an ensemble 
of coupled nonlinear oscillators [22] 


ilfrn ~( J ^n(Pn)4'n ^ni.Pn)4^n + ^ ^ Jnn'^n' 

n' 

+ V D n (Pn)Tn^n » (7) 

where the complex SW amplitudes are defined by 

Ml - iM n v 

tn = , " = . (8) 

V2 M s {M s + M?) 

By writing Eq.(8) as 0^ = y/Pn(t)e^ n ^\ one can see that 
describes the precession of M in the x-y plane, while p n = 

| if/ n | 2 (referred to as the local SW power) is related to the polar 
angle through 6 n = arccos(l - p n ), see Fig. 2. 

The first two terms on the right hand side of Eq.(7) are re¬ 
spectively the nonlinear frequencies a) n (p n ) and damping rates 
T n(Pn) of the n± disk. Both are proportional to the effective 
field y\H ^ ff • z\ acting on each disk. 

For the small precession amplitudes considered here, non¬ 
linear effects are taken into account pertubatively by expand¬ 
ing into powers of p n the frequencies and damping rates, re¬ 
spectively as [22] oj n (p n ) « oj n + vp n , T n (p n ) « aoj n ( 1 + Qp n ). 
Here oj n , v and Q are the coefficients of the expansion of 
(jL> n (p n ) and T n (Pn) to the first order in p n . The third term 
on the right hand side is the interlayer coupling J nn > . In gen¬ 
eral it is a complex quantity, and its phase is related to energy 
gain or dissipation [18]. In the following we will consider the 
simple case of a uniform nearest-neighbour interaction that 
amounts to retain only terms containing if / n± i in Eq.(7) and 
set J nn > = J 6 n ,n±\ . Note also that the imaginary part of J nn > 
and Y n are related since they both stem from the dissipative 
term proportional to a in the LLG equation, Eq. (1). Both 
the nonlinearity and the coupling are due to the dipolar field 
Eq.(3c). The parameters (aj n , v, Q, J nn 0 can be calculated an¬ 
alytically only in some simple cases, but in general they must 
be inferred from micromagnetic simulations [22, 25]. 

We consider the case where the temperature is uniform 
within each disk. In this case, the last term of Eq.(7) de¬ 
scribes thermal fluctuations in terms of the complex Gaus¬ 
sian random variables = ( jf x + ip”) and the nonlinear 
diffusion constant D n (p n ). In the linear regime, the latter 
equals yD , where D is the quantity defined in Eq.(5). In 
the nonlinear regime D n (p n ) depends on T n (p n ), 0 J n (p n ) and 
J nn > and must be fixed consistently to satisfy the fluctuation- 
dissipation theorem. For the single oscillator this amounts to 
fix D n (p n ) oc r n (p n )/aj n (p n )[22]. This is necessary to ensure 
that, for T n = T, the systems approaches a global canonical 
equilibrium at temperature T. Consistently with the small- 
amplitude limit, we assume T n (p n ) = a(L> n (p n ). As a conse¬ 
quence, the noise term in Eq.(7) becomes purely additive and 
Dn(Pn) reduces to a constant D = a in units with Jcb = 1. 


Taking into account the above approximations, we simplify 
Eq.(7) into 

ilpn =(1 + ia) [-y\^n\ 2 ^n ~ - Ji'l'n+l + 

+ y/aTn€ n (i) . (9) 

Upon defining the “Hamiltonian” 

<H = Yj[0J n \^ n \ 2 + V -\iP n \ A + J(r n fn + 1 + « +1 )L (10) 

n 

where (if/ n , hf/* n ) are canonically conjugate variables satisfying 
the Hamilton equations ij/ n = -dfH / di\\j* n , Eq.(9) can be writ¬ 
ten more concisely in the form of a Langevin equation with 
uniform bath coupling a 


i\j/ n = -(1 + ia) 


d<H 

Wn 


+ ^aT n % n . 


In the absence of coupling with the external baths (a = 0), 
Eq.(9) is the well-known DNLS equation [26] which, at vari¬ 
ance with its continuum limit, is not integrable [27]. Such 
equation describes a large class of conservative oscillating 
systems. Some examples include transport in biomolecules, 
Bose-Einstein condensates in optical lattices, mechanical os¬ 
cillators and photonics waveguides [28, 29]. The dependence 
of oj n on the lattice sites introduces an heterogeneity that is 
also connected with the nonlinear version of the Anderson 
tight-binding model [14, 30, 31]. 

In the non dissipative limit, the model admits a further con¬ 
stant of motion besides energy, namely, the total SW power 
P = YjnPn- As a consequence, the thermodynamic equilib¬ 
rium phase-diagram is two-dimensional, and each equilibrium 
state is determined by the energy density 7T/A and the SW 
density P/N [32]. 

In the nonequilibrium regime, the local fluxes of the con¬ 
served quantities are of special interest. By computing the 
time derivatives of the SW power p n and of the local energy 
77^, one obtains the two continuity equations [17, 18, 22, 33, 
34] 

Pn =Jn + Jn + 1 “ Jn > ( lla ) 

U, =J,f + j'L -./,f. (lib) 

where 


Jn =2 ./ n+ \ - if/n)]. 


fn =2/Re 


iW_ 

dilj/n 


OAn +1 - fn) 


(12a) 

(12b) 


are the magnetisation and energy currents associated with the 
hamiltonian coupling between neighboring oscillators. The 
remaining terms 27^ and 27if account for the exchange with 
the reservoirs due to both fluctuations and damping. Their 
steady-state averages are computed using stochastic calculus, 
by evaluating the change of p n and 77„ up to second order in 
the noise term [1]. Asa result, one gets 

(Jn ) = 2aT n - 2 a ([aj„(p n )p n + ./Re (f* n (f„+i + <A«-i)]) 

(13) 
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A similar (but more involved) expression holds for the energy 
fluxes. As a preliminary test, we verified that for a generic 
nonequilibrium stationary state the above definitions of fluxes 
satisfy the local flux balance expressed in Eqs. (11a) and 
(Hb). 

For a system which is driven out-of-equilibrium from its 
boundaries, the local temperature represents a useful observ¬ 
able for the characterisation of the stationary state. Note in 
particular that the quantity T n that appears in Eq.(7) specifies 
the temperature of the phonon bath, which in general does not 
correspond to the temperature of the system [15]. The defini¬ 
tion of temperature for a system of interacting magnetic mo¬ 
ments is not straightforward, since the model Hamiltonian is 
non separable, and one cannot relate temperature to the aver¬ 
age kinetic energy. Within the DNLS formalism, one can use 
the general microcanonical definition of temperature for non 
separable Hamiltonians with two conserved quantities [35]. 
The general expression is nonlocal and rather involved, we re¬ 
fer to Refs. [15, 17, 35] for details. Since we are interested in 
the limit of low temperatures and low amplitudes, we can fol¬ 
low the derivation in Ref. [17] and introduce a simple approxi¬ 
mation of the microcanonical temperature based on a mapping 
of the DNLS equation to a chain of nonlinear coupled rotators 
(XY model). Accordingly, the temperature is approximated 
by 


TxY = fdPn)) ~ (&n) 


(14) 


and it acquires a simple interpretation of the phase- 
fluctuations of the oscillator ij/ n . The function f((p n )) is a 
rescaling factor that depends on the average local power ( p n ). 
One can also show that for ( p n ) 1, f((p n )) = 2 ( p n ). 


IV. MICROMAGNETIC SIMULATIONS 

The micromagnetic simulations were performed with the 
Nmag software [36], using a tetrahedral finite element mesh 
with maximum size of 3 nm, of the order of the Py exchange 
length, and an integration time step is 1 ps. The mesh was 
automatically generated by the Netgen package [37]. 

Each disk of the chain has thickness t - 3 nm, radius R = 
20 nm, and an interlayer distance d = 3 nm. The applied field 
iT ext = IT defines the precession axis of the magnetisation 
along z. The exchange stiffness A = lxlO -11 J/m corresponds 
to that of Py, while the other micromagnetic parameters are 
M s = 0.94 T///o, cc - 8 x 10 -3 and y = 1.873 x 10 11 radxs -1 
xT" 1 . 

The output of the simulations consists of the ensemble of 
the magnetisation vectors {M n (r n , t)}, n = 1,..., 10. Each 
vector depends on the spacial coordinate r n inside the n t h 
disk, and the collective magnetisation dynamics of each disk 
is given by the volume average Eq.( 6 ) [22, 34] 

To have a basis for comparison, we consider first dynam¬ 
ics at zero temperature. Starting with the magnetization uni¬ 
formly tilted 5 with respect to the z axis, the time evolution 
is computed for 30 ns. The SW power spectrum, shown in 



Figure 3: SW spectrum of the system at zero temperature, charac¬ 
terised by 5 dipolar modes. 


Fig. 3, is given by the absolute value of the Fourier transform 
of the collective variable ¥(0 = The spectrum 

consists of five dipolar modes with frequencies (co i,..., 60 ) 5 ) of 
respectively (18,19.7,20.5,21.8,24) GHz. By inspecting the 
spectra of the individual disks, one can see that the mode co\ 
corresponds to the precession of the first and tenth disks, the 
mode C 02 to the second and ninth disks and so on until the fifth 
disk which precesses with frequency see Fig.(l). This in¬ 
dicates that the system has a mirror symmetry around its cen¬ 
ter, due to the fact that each disk behaves as a magnetic dipole. 
Aligning those dipoles in a chain gives a structure where the 
intensity of the dipolar field, which controls the frequencies, 
is symmetric around the center of the chain. 

Let us now discuss the off equilibrium dynamics. We con¬ 
sider the configuration normally adopted to study heat transfer 
in chains of nonlinear oscillators [ 1 ], where the thermal baths 
act only at the boundaries of the system. This setup is differ¬ 
ent from the usual micromagnetic studies of the spin-Seebeck 
effect [ 8 - 11 ], where each spin is coupled to a thermal reser¬ 
voir with a different temperature. The main advantage of our 
choice consists in that it allows to observe the spontaneous 
thermalisation of the system, by probing the XY temperature 
defined in Eq.(14), for the spins that are not directly connected 
to the thermal baths. Experimentally, this could be realised by 
separating the disks with thermal insulating spacers, so that 
energy flows are carried only by the dipolar coupling. 

The simulations at finite temperature were performed start¬ 
ing with the magnetisation aligned along z and evolving the 
system in the presence of the thermal baths for 110 ns. The 
relevant observables were computed after an interval of 70 ns, 
necessary for the system to reach the nonequilibrium steady 
state. The results were time averaged over the last 40 ns and 
then ensemble averaged over 32 samples with different reali¬ 
sation of the thermal noise. 

The high temperature bath T+ ranges between 5 and 30 K, 
while the lower temperature bath T_ is kept fix at 5 K, so that 
AT = T+ - T- is comprised between 0 and 25 K. The local 
currents are always expressed per unit coupling and are thus 
pure numbers. According to our convention, j^ /E refers to 
the current propagating from disk n to disk n + 1 and positive 
(resp. negative) currents propagate from left to right (resp. 
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Figure 4: Profiles of magnetisation (a) and energy (b) currents com¬ 
puted for different values of AT. The profiles are symmetric for 
AT = 0, and then become strongly asymmetric. The lines are guide 
to the eye. 


from right to left). 

In contrast to previous studies of the off-equilibrium DNLS, 
here the dissipation is both at the edges and in the bulk of the 
system. As a consequence, the currents do not have a flat 
profile in the bulk, but decrease exponentially along the chain. 
In the thermodynamic limit, this system is thus an insulator. 

Figs. 4 (a) and (b) show respectively the profiles of magneti¬ 
sation and energy currents, for different values of AT. The two 
currents have similar profiles, and they do not vanish when 
AT = 0 (black dots). In both cases, they are symmetric with 
respect to the zero-current axis. This behaviour is due to the 
fact that the currents generated by the two baths travel in op¬ 
posite directions and decrease because of the damping, van- 
ishng in the middle of the chain. Note that, although the local 
currents do not vanish at equilibrium, the total current is zero. 



Figure 5: Profiles of the magnetisation currents for 77 = 0 K and dif¬ 
ferent values of T+. The right-going current decreases exponentially 
along the chain. The dashed lines are fit with Ae ~ x/x °, for different 
values of the parameters (A,x 0 ). The energy currents have similar 
profiles (not reported). 


1 - 1 - i - 1 - i - 1 - i - 1 - 

u 0 6 12 18 24 

AT (K) 

Figure 6: Currents j M/E = jf /E + jg /E injected from the reservoirs vs 
the temperature difference AT. The currents grow linearly with AT 
in the low temperature regime, and then reach a plateau. The lines 
are guides to the eye. 


Increasing T+ leads to the increase of the right-going currents, 
and moves the point where the local currents vanish. 

In all these cases, the system behaves essentially as sink: 
the currents injected from the baths are dissipated in the bulk, 
so that there is no net transport. A situation where trans¬ 
port occurs can be observed in Fig. (5), which shows the case 
where 77 = 0 and (the right-going) current is injected only by 
T +. The current decreases exponentially, and for high enough 
value of AT remains positive until the end of the system. Note 
that in the thermodynamic limit this system is an insulator, due 
to the damping in the bulk. 

Note that the local currents do not grow indefinitely with 
AT, but they saturates around AT « 20 K. This can be seen 
more clearly in Fig. 6, where the currents injected from the 
reservoirs j M/E = j^ /E + j^ /E grow linearly with AT until 
they reach a plateau for AT « 20 K. 

This phenomenon originates from the fact that the LLG 
equation is nonlinear, and consequently the frequency spec¬ 
trum of the system is temperature dependent. Fig. 7 shows 



Figure 7: Finite temperature power spectra of disks 1 and 2. (a) 
At low temperatures (T± = 5 K) regime the peaks overlap and the 
currents are proportional to AT. (b) At high temperature (T+ = 30 
K) the spectrum of disk 1 shifts towards high frequency, reducing its 
overlap with disk 2. In this condition, the current does not increase 
anymore with AT. 
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Figure 8: (a) Profiles of SW powers p n for different values of AT. 
(b) Profiles of spin temperature T XY for different values of AT. The 
inset displays a magnification of the profiles between the disks 2 and 
8, which show the spontaneous thermalisation of the chain. The lines 
are guide to the eye. 


the power spectra of the disks 1 and 2. Panel (a) displays the 
low temperature regime, with T+ = T- = 5 K. Thermal fluc¬ 
tuations excite all the modes of the system, and both disks 
display a broad peak between 17 and 25 GHz. The current 
increases linearly with AT as far those peaks overlap. Pan¬ 
els (b) show the high temperature case, where disk 1, directly 
connected to the hot bath, increases its frequency and does 
not overlap with disk 2. This situation is a manifestation of 
stochastic phase synchronisation (that is, the control of syn¬ 
chronisation through temperature) in the propagation of heat 
current. A similar mechanism is at the basis of spin and ther¬ 
mal rectifiers [33, 34, 38^10]. 

The profile of the local SW powers p n is displayed in Fig. 8 
(a). The powers reach the maximum at the edges, directly con¬ 
nected to the baths, and decay along the chain. As expected, 
their profile is symmetric for AT = 0, and becomes strongly 
asymmetric when the temperature difference is finite. Note 
that the powers increase with with temperature until AT ^15 
K and then remain roughly constant. In this high temperature 
regime, p\ keeps increasing because of thermal fluctuations, 
but p 2 decreases. This is due to the fact that energy remains 
confined in the first disk due to the de-synchronization, as pre¬ 
viously discussed. 

Fig. 8(b) shows the spin temperature profiles T X y of the 
chain, computed for different values of AT. The temperature 
is higher at the boundaries, where phase fluctuations are given 
by the direct contact with the baths, and then decreases dra¬ 
matically towards the center of the chain, as can be seen from 
the inset. 



Figure 9: Simulation of Eq. (9) with J = 0.1, a = 0.008, v = 1, 
(a)\, 0 ) 2 , 0 )$, a) 4 , (Os) = (1,1.09,1.15,1.21,1.33) and different applied 
temperature differences AT . The temperature of the leftmost reser¬ 
voir spaces in the interval T+ = [0.5, 1.25], while the rightmost one 
is kept fixed at T_ - 0.5. 


V. COMPARISON WITH THE OSCILLATOR MODEL 

The micromagnetic simulations indicate that the system has 
some form of mirror symmetry around its center. Thus we 
model this by choosing the linear frequencies in model (9) 
such that a)® = 0J° N _ n . A relatively small coupling / = 0.1 has 
been employed. As in Section IV, the fluctuating forces are 
applied only at the boundaries of the system, i.e. we choose 
T n = \T+S n , i + T- S ntN ]; 



Figure 10: Simulation of Eq. (9), SW profiles and microcanonical 
temperature 7^ for different applied temperature differences AT. The 
parameters of the simulation are the same as in Fig. 9. For each 
lattice site n, the local temperature 7^ is calculated by evaluating the 
general definition [15,17] on a triplet of sites centered around n. The 
inset of the lower panel shows a comparison between T^ and T X y for 
the profile with AT = 0. 
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Figure 11: Simulation of Eq. (9) Boundary magnetisation (a) and 
energy (b) currents toward the two left and right reservoirs versus the 
temperature difference AT for T_ - 0.05, other parameters given in 
the text. Panel (c) shows the corresponding excess fluxes calculated 
with Eq.(15). 


As a test for stationarity we evaluated the average currents 
from the heat baths 2^(27^) an d > are indeed vanish¬ 

ing within statistical accuracy. 

In Figs. 9 and 10 we show the current and SW power pro¬ 
files along the chain. Comparing with the corresponding fig¬ 
ures 4(a,b) and 8(a), the qualitative agreement is good. The 
lower panel of Fig. 10 reports the temperature profiles 7^ cal¬ 
culated with the exact microcanonical expression defined in 
Refs. [15, 17]. Such profiles are in good agreement with the 
phase temperature T X y (see the inset) and qualitatively repro¬ 
duce the temperature profiles calculated within the micromag- 
netic framework in Fig. 2. Altogether, this confirms that defi¬ 
nition (14) is a sensible approximation in the present setup. 

Finally, in panels (a) and (b) of Fig. 11 we show the sta¬ 
tionary boundary fluxes and 27jf • Nicely, they reproduce 
the saturation effect observed in Fig. 4 with the micromag- 
netic simulations. It should be noted that the currents do not 
vanish for AT = 0: this is simply because the system is off- 
equilibrium also in this case due to the presence of the dissi¬ 
pation in the bulk of the chain. For the very same reason, the 
two currents flowing at the boundary are not equal. In order to 
single out the purely transport contribution from the dissipa¬ 
tive one in the boundary fluxes, we define an excess boundary 
flux as 


® M ’ £ =jM’ e (AT) -J^ e (0) , (15) 

that takes into account the amount of magnetisation/energy 
which is transported per unit time with respect to the purely 
dissipative symmetric profile at AT = 0. The behaviour of 
<& me (AT) is shown in Fig. 11(c). Both the excess fluxes have 
positive sign, meaning that there is a net transport of energy 


and magnetisation from left to right. Moreover they show the 
same saturation effect observed in panels (a) and (b). 


VI. CONCLUSIONS 

In this work we have illustrated the spin-Seebeck effect in 
a small array of coupled magnetic nanodisks and we have 
related to general transport properties of out-of-equilibrium 
chains of nonlinear oscillators. The dynamics of the disk 
chain has been investigated by means of micromagnetic simu¬ 
lations and compared with the effective model, Eq. (9). There 
is indeed a good qualitative agreement between the two ap¬ 
proaches, and a good physical insight can be achieved from 
the coupled oscillator model. 

Another remarkable result, is that the relationship with the 
XY model provides a simple prescription for computing the 
local phase temperature in the micromagnetic simulations via 
Eq. (14). Such observable plays a relevant role in our setup, 
since it allows to quantify thermal fluctuations inside the sys¬ 
tem, i.e. for macrospins that are not directly connected to 
the external reservoirs. This last issue is of major importance 
for non-standard Hamiltonians like the DNLS one, where ki¬ 
netic and potential energies are not separated. Indeed, the 
XY approximation allows introducing the simple kinetic ex¬ 
pression T X y for the temperature, that can safely approximate 
the microcanonical one 7^. This is of practical importance, 
considering that the microscopic definitions of T and fi are 
pretty much involved for a non separable Hamiltonian, like 
the DNLS one. 

The presence of two coupled currents is related to the ex¬ 
istence of two thermodynamic forces [41^13] (AT, A//), the 
latter being the difference of chemical potential [15, 17]. In 
this work we limited to the case where Afi = 0. Actu¬ 
ally, chemical potential gradients can be easily accounted 
for within the DNLS language, by adding terms of the form 
ia/i n i// n to Eq.9 [17]. For our system of precessing spins, 
this is interpreted as a torque that compensates the damping 
and controls the magnon relaxation time towards the reser¬ 
voirs [15, 17, 34]. This can be experimentally realised trough 
spin transfer torque [22, 44, 45]. In this case the spin-torque 
induces a nonequilibrium dynamics and may lead to self- 
sustained oscillations, thus opening a new realm of transport 
phenomena. We plan to investigate those setups in the next 
future. 
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